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ABSTRACT 

A program for molecular calculations with B functions is reported and its performance is 
analyzed. All the one- and two-center integrals, and the three-center nuclear attraction 
integrals are computed by direct procedures, using previously developed algorithms. The 
three- and four-center electron repulsion integrals are computed by means of Gaussian 
expansions of the B functions. A new procedure for obtaining these expansions is also 
reported. Some results on full molecular calculations are included to show the capabilities 
of the program and the quality of the B functions to represent the electronic functions in 
molecules. 

1 Introduction 

Although there is a general agreement in the best quality of exponential type orbitals 
(ETO) for describing the electronic function in atoms and molecules the difficulties 
in solving some molecular integrals (in particular the three- and four-center two-electron 
repulsion integrals) have strongly restricted their use in molecular calculations, in benefit 
of the more easy to handle Gaussian type orbitals (GTO). 

In a previous article 0, a program for molecular calculations with Slater type orbitals 
(STO) has been reported. In that program, some of the molecular integrals (all the 
one- and two-center and the three-center nuclear attraction) were computed by direct 
algorithms, and the remaining ones by means of Gaussian expansions of the STO. The 
aim of that program was to provide a bridge and a reference for further works in this 
field, mainly to those dealing with the direct calculation of all the integrals. The present 
work has been developed in the same spirit, using in this case a different type of ETO, 
the reduced Bessel functions (B functions) that have been already proposed by Steinborn 
et al. m |[] for molecular calculations. 
In the real case, the B functions are defined as: 

B™ t (a, r) = [T +l (n + l)!]" Vi/a(«r) (ar) 1 z^O, </>) (1) 
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where 

k r (z) = y| z* K 7 (z) (2) 

are the corresponding reduced Bessel functions, with being the Macdonald function 
H, zf 1 are real spherical harmonics: 



(2/ + l)(/-|m|)! 
\\ (1 + S mi0 ) 7T (l+\m\)\ 



-l) m pl m| (cos^)$ m (0) (3) 



P ; ' m ' being the corresponding Legendre functions |J, and the functions $ are defined as: 

^ / ,\ _ | cosm^ for m > 

— j sin | m |0 for m < ^ ' 

The simplicity of the Fourier transform of the B functions has been proposed as an advan- 
tage to be exploited in the development of algorithms for molecular integrals calculations 
with this kind of ETO || |. 

In the next sections, the algorithms used in the program will be summarized and some 
results will be reported to show the quality of the B functions for molecular calculations 
and to compare with both the commonly used Gaussian functions and with the STO. 



2 Algorithms 

The present program for molecular integrals with B functions closely follows the scheme 
of that previously reported for the STO ||. The one- and two-center and the three- 
center nuclear attraction integrals are computed by direct procedures previously developed 
|| f|, [7], H, D, [L(J and the remaining ones are computed by means of Gaussian expansions 
of the B functions. 

The integrals computed with direct procedures are obtained with an average cost about 
0.1 milliseconds per integral actually computed measured on a Digital AlphaServer 8400, 
for an accuracy of at least twelve decimal places. 

The remaining integrals have been computed by using Gaussian expansions, as mentioned 
above. The expansions have been obtained by a procedure different than the usual least- 
squares method. This procedure, that is summarized in the appendix, is much easier to 
apply than the least-squares, and yields expansions that tend to better approximate the 
tails of the functions. 

The integrals between the Gaussian primitives are evaluated with the algorithm proposed 
by Saunders [□], slightly modified, and a test for avoiding calculation of negligible 



Gaussian contributions to a given set of integrals is also included. 
The program for the calculation of the integrals has been implemented in a modular way, 
so that the subroutines corresponding to each type of integrals can be easily replaced by 
others when required. The main structure can therefore be kept unaltered when trying 
new algorithms for a given type of integrals. 

The package includes routines for the direct minimization of the energy in RHF and 



ORHF calculations [13[ but, since the integrals are stored in external files in a very easy- 
to-handle way, the interface with other standard programs for minimizing the energy 
(both in Hartree-Fock and post-Hartree-Fock calculations) is straightforward. 
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3 Results 



We have first analyzed the quality of the Gaussian expansions for representing the B 
functions. In Table 1 we collect the least-squares error of the expansions for different 
principal quantum numbers in the functions. We want to recall here that the least-squares 
error has not been used as criterion for attaining the present expansions. 
As it can be seen in the table, the quality of the expansions improves as the n quantum 
number increases. This is so because the Gaussian expansions obtained with the current 
procedure tend to better reproduce the tail of the functions, but present some problems 
to reproduce the peak in the origin for the functions with n — 1. This causes a loss of 
accuracy in the integrals involving these latter functions. In consequence, the program 
uses expansions optimized by minimizing the least-squares error for k\/2 (which are the 
same as those of the STO Is £L4|), and expansions obtained with the current procedure 
for the rest. 

Then, we have carried out full RHF calculations on several systems to test the performance 
of the program and to analyze the qualities of the B functions for molecular calculations. 
The corresponding exponents of the single-zeta and double-zeta basis sets used in the 
calculations are collected in Tables 2 and 3. These basis sets are optimized for atoms and 
have been taken from references Jl5, 16]. The geometries are summarized in Table 4. 
In Tables 5 and 6, the energy values and the computational cost are analyzed for ex- 
pansions of different lenght with both types of basis sets. The computational cost of the 
integrals calculated by direct procedures is included in the second column of each table, 
and columns 3 to 8 collect both the value of the electronic energy and the time required 
for the integrals computed with the Gaussian expansions for the B-10G, B-15G and B- 
20G expansions. As it can be seen in these tables, the direct procedures enable us to 
obtain highly accurate integrals in a very fast way. This should spur for a further search 
on direct procedures for all the integrals. 

Nevertheless, an accuracy sufficient for testing purposes can be attained with a moderate 
cost by using not very long expansions (such as the B-10G) of the B functions. In fact, 
the results obtained with these expansions are also sufficiently accurate to give an idea 
about the quality of the B functions in molecular calculations. In Table 7 a comparison 
between the results obtained with B and STO basis sets of same length is made. The 
table clearly illustrates that, despite that the energy values with STO are slightly lower 
in the single zeta basis sets, the double zeta of both types yield results of similar quality, 
thus confirming the capabilities of the B functions for high-quality molecular calculations. 



4 Conclusions 

A program for computing molecular integrals with B functions has been implemented 
and tested. The program combines techniques of direct computation for the one- and 
two-center integrals and the three-center nuclear attraction integrals with others based 
on Gaussian expansions for the remaining ones. Direct procedures yield highly accurate 
results at low computational cost. The Gaussian expansions lead to a lower ratio accu- 
racy/cost, but still acceptable for calculations in medium size systems. The expansions 
obtained by the moments procedure reported herein are good enough for all cases except 
the k\/2 function because of the peak in the origin in this latter. An expansion based in 
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the least-squares procedure is preferable in this case, and has been used in the program. 
Finally, basis sets of B functions have a good quality for reproducing the electronic wave 
function in molecules, comparable with that of STO. 



5 Appendix: Gaussian expansions of B functions 

The procedure used to obtain Gaussian expansions of the B functions is based on the 
fitting of the value of the functions at the origin and the moments of the unnormalized 
spherical part of the B-function, i.e., an expansion: 



N 



kj. 



-i/2(r. 



~ N ,Ci e 



(5) 



is chosen with the requirement that a given number of the moments of the reduced Bessel 
functions are reproduced: 



/*oo A /*oo 

/ dr r fc_1 k n _ 1/2 (r) = V c { / dr 
Jo a Jo 



r k-l e -a t r* k=l,2,... 



(6) 



Replacing the integrals by their final values [IS, 20 1, the following equalities are obtained: 
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(7) 



Introducing the definitions: 



Qo = fcn-l/2(0) 
Q k = 2 fc+ - 1 / 2 r ^ + n) 



and 

Eq. (0) can be rewritten as: 



n = 



N 



Qk = T, c i r " k = 0,l,2N-l 



i=l 



(9) 
(10) 

(11) 



Since an expansion of length N implies 2N unknowns, Cj, r iy a total of 2N moments 
(k = 0, 1, ...2N — 1) can be exactly reproduced. This can be accomplished by solving a 
system with 2N equations that can be written in matrix form: 



/ Qo \ 

Qi 

Qn-i 
Qn 

V Q2N-1 j 
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(12) 
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It should be noted that in this system, both the q and the are unknown. To solve it, 
the system must be partitioned in two, one corresponding to the upper half equalities and 
the other to the lower half ones. By defining: 



Q 



/ Qo \ 
Qi 

V Qn-i J 



Q' 



/ Qn \ 

Qn+i 

V Q2N-1 J 



( c l \ 

c 2 

V C N J 



(13) 



/ 1 

7*1 



1 

r-2 
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( n 



N 
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N-l 



T 



it is clear that Eq. (|T2| ) is equivalent to: 

Q = Tc 

Solving for the first one: 



Q' 

T 1 Q 



'1 



T c 



'2 
r N+l 
<2 



N-l 



' N 

N 



,iV-l 
N 



(14) 



and replacing in the second one, it follows: 

Q' = T' T 1 Q 



(15) 
(16) 

(17) 



To attain directly the general solution of this latter equation can be rather difficult. 
However, one can start by the first cases (TV = 1,2,...) and proceed by induction. In this 
way, it can be found that the set {rj}^ that fulfills the equation conincides with the 
roots of the iV-th degree polynomial: 



N 



i=0 



(n) 







AN) 



where the coefficients, are the solutions of the linear system: 
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Q2N-2 j 

IH), the q can be readily attained by 



Finally, once the are obtained by solving Eq. 
Eq. (0) and the values of ctj follow from (p!0|). 

Gaussian expansions of the sigma part of the B functions ranging from a single Gaussian 
to 20 Gaussians have been obtained. The standard limits for the three- and four-center 
integrals are N — L < 9. 
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Table 1: Least-squares errors , A 2 , in the Gaussian expansions of k n ^i/ 2 (r) obtained with the method 
of moments. 



Number of 










n index of the k n -i/2 function 






Gaussians 
























1 


2 


3 


4 


5 


6 


7 


8 


1 


2.7( 


-1) 


5.9(-2) 


2.5( 


-2) 


1.4(-2) 


8.9(-3) 


6.1(-3) 


4.5(-3) 


3.4(-3) 


2 


3.2( 


-2) 


3.6(-3) 


8.7( 


-4) 


3.1(-4) 


1.3(-4) 


7.0(-5) 


3.9(-5) 


2.4(-5) 


3 


7.1( 


-3) 


4.6(-4) 


7.2( 


-5) 


1.8(-5) 


5.6(-6) 


2.2(-6) 


9.5(-7) 


4.6(-7) 


4 


2.0( 


-3) 


8.4(-5) 


8.8( 


-6) 


1.5(-6) 


3.7(-7) 


1.1 (-7) 


3.9(-8) 


1.5(-8) 


5 


7.0( 


-4) 


1.9(-5) 


1.4( 


-6) 


1.8(-7) 


3.4(-8) 


8.0(-9) 


2.3(-9) 


7.4(-10) 


6 


2.8( 


-4) 


5.0(-6) 


2.7( 


-7) 


2.7(-8) 


3.9(-9) 


7.3(-10) 


1.7(-10) 


4.7(-ll) 


7 


1.2( 


-4) 


1.5(-6) 


6.0( 


-8) 


4.6(-9) 


5.3(-10) 


8.2(-ll) 


1.6(-11) 


3.7(-12) 


8 


5.7( 


-5) 


4.9(-7) 


1.5( 


-8) 


9.0(-10) 


8.4(-ll) 


l.l(-ll) 


1.8(-12) 


3.5(-13) 


9 


2.8( 


-5) 


1.8(-7) 


4.1( 


-9) 


2.0(-10) 


1.5(-11) 


1.6(-12) 


2.2(-13) 


3.8(-14) 


10 


1.5( 


-5) 


6.8(-8) 


1.2( 


-9) 


4.8(-ll) 


3.0(-12) 


2.7(-13) 


3.2(-14) 


4.6(-15) 


11 


8.3( 


-6) 


2.8(-8) 


4.0( 


-10) 


1.2(-11) 


6.5(-13) 


5.0(-14) 


5.0(-15) 


6.4(-16) 


12 


4.8( 


-6) 


1.2(-8) 


1.4( 


-10) 


3.5(-12) 


1.5(-13) 


9.9(-15) 


8.7(-16) 


9.6(-17) 


13 


2.9( 


-6) 


5.4(-9) 


4.9( 


-11) 


1.0(-12) 


3.8(-14) 


2.1(-15) 


1.6(-16) 


1.6(-17) 


14 


1.8( 


-6) 


2.6(-9) 


1.9( 


-11) 


3.3(-13) 


1.0(-14) 


4.9(-16) 


3.3(-17) 


2.8(-18) 


15 


1.1( 


-6) 


1.3(-9) 


7.4( 


-12) 


1.1(-13) 


2.9(-15) 


1.2(-16) 


7.0(-18) 


5.3(-19) 


16 


7.3( 


-7) 


6.3(-10) 


3.0( 


-12) 


3.7(-14) 


8.6(-16) 


3.1(-17) 


1.6(-18) 


1.1(-19) 


17 


4.9( 


-7) 


3.3(-10) 


1.3( 


-12) 


1.4(-14) 


2.7(-16) 


8.5(-18) 


3.8(-19) 


2.3(-20) 


18 


3.3( 


-7) 


1.8(-10) 


5.8( 


-13) 


5.1(-15) 


8.7(-17) 


2.4(-18) 


9.6(-20) 


5.1(-21) 


19 


2.3( 


-7) 


9.8(-ll) 


2.6( 


-13) 


2.0(-15) 


2.9(-17) 


7.1(-19) 


2.5(-20) 


1.2(-21) 


20 


1.6( 


-7) 


5.5(-ll) 


1.2( 


-13) 


7.9(-16) 


1.0(-17) 


2.2(-19) 


6.9(-21) 


3.0(-22) 
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Table 2: Singlc-Zcta BTO exponents. 



Orbital 


Zn 


s 


B 


C 


N 





F 


H 


Is 


28.979194 


15.396775 


4.649767 


5.636105 


6.621925 


7.607778 


8.593356 


1.000000 


2s 


9.212368 


4.468108 


1.076139 


1.346562 


1.612481 


1.885508 


2.154463 




3s 


4.615722 


1.723750 


1.226030 


1.581274 


1.929475 


2.238550 


2.561510 




4s 


0.966290 


5.987867 














2p 


13.015418 


1.684294 














3p 


4.754359 
















3d 


4.660219 

















Table 3: Doublc-Zcta BTO exponents. 



Orbital 


Zn 


S 


B 


C 


N 


O 


F 


H 


Is 


41.443334 


16.517740 


7.809120 


9.121695 


7.308276 


8.306577 


9.356483 


1.200000 


Is 


29.336793 


11.035952 


3.907381 


5.156559 


4.211759 


5.330563 


6.202314 


1.000000 


2s 


23.592780 


5.492259 


4.173803 


5.447762 


3.434959 


2.910067 


3.089333 




2s 


12.264347 


4.467502 


1.204488 


1.495602 


1.724930 


1.873454 


2.040237 




3s 


11.550144 


2.255074 


2.213227 


2.726540 


3.240492 


3.686133 


4.174972 




3s 


5.617435 


1.566569 


1.004047 


1.255198 


1.496401 


1.655654 


1.847191 




4s 


4.425408 


9.438960 














4s 


1.140164 


5.089873 














2p 


19.000918 


2.008415 














2p 


11.734571 


1.206954 














3p 


5.624098 
















3p 


3.506195 
















3d 


7.339257 
















3d 


3.134250 
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Tabic 4: Geometries used for molecular calculations 



Molecule 


Geometry 


Bond distances and angles 


BH-) 


Planar 


1 > 


O or: 

z.zo 






See fig. 1 




2.ZOU13 


tx BH — z.ooUo 1 






1 > 


o.oo4oU 


iffBB = iiy 


CH 4 


RpEnil^r totr^borlron 


Kcjj = 


z.0665 




— £ £ 


Linear 


1 > 

RffC = 


z.UUz 


K.CC = ^-zol 


C2H4 


Planar 


Rffc = 


2.02203 


Rcc = 2.55116 








= 120° 




C2H6 


Alternate 


R-ffc = 


2.08250 


Rcc = 2.91588 






^HCH = 


= 109.32° 


Lhcc = 109.62° 


HCN 


Linear 


R-ffc = 


2.0 


R CN = 2.187 


SF 6 


Regular octahedron 


Rsf = 


2.88769 




Zn 3 


Equilateral triangle 


R-ZnZn 


= 5.03593 





11 Distances in a.u. and angles in degrees. 



Table 5: Electronic Energy RHF with Singlc-Zeta basis set of BTO. a 



Molecule 


time(direct.alg.) 


10G 


time 


15G 


time 


20G 


time 


BH 3 


0.04 


-26.2978760 


1.6 


-26.2978758 


6.1 


-26.2978762597 


16.4 


B2H 6 


0.6 


-52.591815 


41.3 


-52.591817 


156.9 


-52.5918199842 


410.4 


H 2 


0.1 


-75.6136703 


0.5 


-75.61366972 


1.9 


-75.6136697787 


5.0 


HCN 


0.2 


-92.4898386 


2.6 


-92.48983974 


9.7 


-92.4898398408 


25.0 


CH 4 


0.04 


-40.0629604 


3.5 


-40.0629605 


13.7 


-40.0629611447 


36.3 


C2H2 


0.3 


-76.5372092 


6.3 


-76.5372103 


24.4 


-76.5372106538 


61.6 


C2H4 


0.4 


-77.7418899 


19.8 


-77.7418903 


74.8 


-77.7418909888 


192.4 


C2H6 


0.6 


-78.970476 


40.4 


-78.9704946 


155.8 


-78.9704966047 


404.1 


SF 6 


9.9 


-988.79773 


598.7 


-988.79379 


2277.8 


-988.7938448335 


5836.4 


Zn 3 


82.4 


-5306.54329 


275.4 


-5306.54320060 


614.4 


-5306.5432005659 


1679.2 



11 Energy in Hartrees. 
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Table 6: Electronic Energy RHF with Doublc-Zeta basis set of BTO. 



Molecule 


time(direct.alg.) 


10G 


time 


15G 


time 


20G 


time 


BH 3 


0.5 


-26.373450 


23.1 


-26.3734463 


90.0 


-26.3734465526 


234.7 


B2Hg 


5.3 


-52.757514 


595.0 


-52.757528 


2290.9 


-52.7575305959 


5967 .9 


H 2 


0.3 


-75.998708 


7.3 


-75.99870643 


28.0 


-75.9987064604 


73.0 




9 k 




oD.4 


09 80.1 1 ARA 


1 ^ fi 

loD.O 


09 801 1 A SRSRS 

-yz.oui 1400000 


oOU.U 


CH 4 


0.9 


-40.182804 


51.6 


-40.1828002 


202.8 


-40.1828008395 


532.5 


C2H2 


3.0 


-76.773130 


89.1 


-76.7731380 


336.4 


-76.7731382182 


872.2 


C2H4 


4.2 


-77.99422 


281.3 


-77.9942558 


1070.7 


-77.9942574766 


2791.2 


C2H6 


5.1 


-79.197399 


589.0 


-79.1974595 


2255.7 


-79.1974630394 


5853.2 


SF 6 


141.5 


-993.7002 


8683.5 


-993.70374 


33498.9 


-993.7038527430 


87495.1 


Zn 3 


744.6 


-5332.8974 


2007.0 


-5332.89797 


6802.8 


-5332.8979948489 


16198.3 



Energy in Hartrees. 



Table 7: RHF electronic energy with both Slater and B basis sets. a 



Bessel(SZ) @ Slater(SZ) @ Besscl(DZ) M Slatcr(DZ) p 



BH 3 

B2HQ 

H 2 

HCN 

CHi 

C2H2 

C2H4 

C2HQ 

SF 6 

Zn 3 



-26.297876 
-52.591815 
-75.613670 
-92.489839 
-40.062960 
-76.537209 
-77.741890 
-78.970476 
-988.79773 
-5306.54329 



-26.319008 
-52.630434 
-75.687535 
-92.572795 
-40.101608 
-76.596599 
-77.810378 
-79.041644 
-989.887733 
-5313.330879 



-26.373450 
-52.757514 
-75.998708 
-92.801148 
-40.182800 
-76.773130 

-77.99422 
-79.197399 

-993.7002 
-5332.8974 



-26.375009 
-52.758232 
-76.000535 
-92.833977 
-40.183596 
-76.803955 
-78.004891 
-79.199828 
-993.705704 
-5332.907860 



a Energy in Hartrees. 
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